叶绿体/线粒体在线注释网站GeSeq使用简介
上篇介绍了一个线粒体在线注释网站,MITOS,今天呢白鱼小编继续给大家推荐另一个在线注释网站,GeSeq。该网站即可用于注释线粒体,也可用于叶绿体。GeSeq主要以同源注释的方式,快速准确地实现细胞器基因组的注释,特别是叶绿体基因组。
GeSeq使用简介
GeSeq界面功能
打开GeSeq网站链接:https://chlorobox.mpimp-golm.mpg.de/geseq.html,上传序列,调节参数。以下是对GeSeq界面的基本简介。
左侧上传序列(fasta文件),并选择基本选项,包括序列类型、注释类型并设置参数信息等。
GeSeq主要依靠同源注释,将你上传的基因组序列和参考基因组序列作比对,来确定你的基因组中的结构区。同源比对的方法有blast和hmmer两种,参考基因组很近缘的情况下,blast足够使用了。
blast和hmmer既能注释蛋白编码区(cds区),也能注释非编码RNA区(tRNA、rRNA)。蛋白编码cds区目前只能通过同源注释获得。对于非编码RNA区,除了同源比对,GeSeq还提供了3种从头预测非编码RNA的方法,其实就是调用3种非编码RNA预测软件来实现。如果同源比对未能找到全部的非编码RNA,那么就可以通过这3种方法来尝试获取。
正常情况下,先将这3种方法关闭,若同源注释缺少结果时再行开启也不迟。有一点需注意,这3种方法之间相互独立,和同源注释方法也独立,若开启多种功能时,就会导致出现重复的RNA区域(因为各软件各自独立运行获得结果,GeSeq自己并不会去重),这时候还需手动在结果文件中删除重复、冗余的结果。
中间部分,导入参考基因组序列(同种或近缘物种的叶绿体/线粒体基因组序列)。若想保证GeSeq正常运行,这里必需要载入至少一个参考基因组。
GeSeq提供了在线导入和本地导入两种方式。无论使用哪种,参考基因组肯定是越近缘越好,这样就能保证准确性和完整性;否则结果中将会出现很多基因无法准确注释,主要体现在边界区无法准确确定,覆盖度不完全,长度过短,未得到有效的编码氨基酸序列等(这时候你还需要手动去鉴定,会很麻烦,这也是很多测序公司不愿意接收缺乏近缘参考基因组的细胞器的注释项目的原因,注释繁琐是一方面,有时候可能连基因组组装都是问题)。GeSeq和NCBI的数据库相连的,所以一般来讲我们输入的参考基因组的NCBI登记号能够很轻易找到并在线导入。
前面的选项设置好之后,若无问题,右侧点击提交,等待一会儿即可出结果。注释结果会在服务器中保存一段时间(时间不是长达几周,或者你没有手动删除的话),期间若你浏览器关闭后再次打开GeSeq网站,会自动加载来的结果。
序列上传注释及结果说明
下面我们以某植物叶绿体基因组为例,演示该网站的使用。该基因组以fasta格式存储,文件名称“test.fasta”,里面包含一条完整的环状叶绿体序列,序列id为contig1。(我这儿就不把这段序列上传了,后面的步骤,大家使用自己的序列模仿着来即可;若是觉得叶绿体注释繁琐,可以换个动物线粒体尝试)
对于该叶绿体序列,我们的参数设置如下。
最后提交,等待出结果。
点击查看文件,在打开的文件中点击“Download”可将结果下载到本地。
Genbank注释结果文件。
GFF注释结果文件。
基因核酸序列fasta文件。
检查结果的完整性,以及后续调整
是不是觉得注释完成了?不,我们还需检查结果的完整性、可靠性。
GeSeq的结果文件需要后续调整后才能使用,不信你运行一个后点开gff文件查看,里面的顺序乱糟糟的,你得调整……不过这倒不是主要的问题,首先要查看的就是我们的基因组中,蛋白编码基因和非编码RNA区域是否有遗漏的,或者没注释完全的。
就我个人而言,我习惯下载下来Genbank注释文件后,在Genbank文件中修改。因为Genbank文件中的内容全面,而且改好这一个后可以直接基于gbk文件再转成gff文件、基因核酸序列fasta文件、基因编码氨基酸序列fasta文件、上传NCBI用的tbl文件等。这样也就相当于对其它的文件同时作修改了,就无须再单独修改原始的gff、fasta文件之类的了,快捷而高效。
我继续举一些经验例子吧。我习惯的做法,打开参考基因组的gbk文件,以及GeSeq注释得到的gbk、文件,两者之间相互比较(眼睛看,嗯嗯,一定要仔细啊),看自己的注释结果和参考基因组的在哪里出现了不同。
备注:为了更方便比较,推荐在基因组注释前,将你的叶绿体基因组序列的起点调整为与参考基因组序列一致。(调整起点只对环状序列适用)
一般来讲,近源物种的叶绿体基因组还是非常保守的(无论种类还是碱基组成),所以大部分基因、RNA等位置都是一致的。对于不一致的区域,根据注释得到的基因(或RNA)名称,看参考基因组中这段位置是否是这个基因,以及上下游的基因是否一致。如果这儿缺少基因了,不妨重新调整参数注释(比方说blast参数调宽松一些,开启hmmer比对,若有未得到的RNA,还可启用从头预测的方法);不行的话还需结合其它的细胞器注释工具获得(GeSeq网站中也列出了一系列可用的替代工具,https://chlorobox.mpimp-golm.mpg.de/Alternative-Tools.html),或者手动查找,比方说本地blast等,虽然繁琐但也是没办法的事情。若出现多出来的基因,先看下这个基因的blast结果是否严谨,如果coverage也就40%的那种,那么基本就可断定这条注释是错误的,需要剔除。
经常会出现相差一些碱基的情形,对于非编码RNA还好说(这个只需找到核酸序列就行了),相差几个碱基也能说得过去;关键在于蛋白编码基因区域的确定(这个还需给出完整的编码氨基酸序列,需要明确起始、终止密码子的位置),这个就严格了,相差几个碱基就得不到有效的编码产物。GeSeq注释结果中,经常会遇到这样的情形。参考基因组与你的基因组亲缘关系越远,就越容易出现未能注释完全的基因结构,如下示例就是其中情形之一。
那么,碰到这种情况时,怎么手动补全呢?我继续分享经验。
这个时候,我们就要根据参考基因组和自己基因组中这个基因的位置,手动推测基因边界,并作鉴定。就以上图为例吧,参考基因组中,该基因“rpl16”的cds区所在区域为“join(182497..182505,183506..183904)”,长度为“9bp+399bp”;我们的基因组中,得到的“rpl16”的cds区位于“183624..184022”,长度为“399bp”。这时候我们就很容易推测,GeSeq注释时只注释到了后面那一段cds区域,未能有效识别前一段cds以及内含子结构,很可能由于这段cds区长度过短所致。由于叶绿体基因组非常保守(至少近缘物种的是这样,变异很小),这就很有利于我们准确推测前面那小段cds区的位置:对于后面那一段cds区,参考基因组中的位置为“183506..183904”,我们的基因组中的位置为“183624..184022”,起始和终止位置“相差118bp”,即可知,前面小段cds区的“相差间隔”可能也是118bp,根据参考基因组中的位置“182497..182505”,推测其在我们的基因组中的位置“182615..182623”。
好了,这只是推测,怎么作验证呢?最简单的方法就是把参考基因组中“join(182497..182505,183506..183904)”,以及我们的基因组中“join(182615..182623, 183624..184022)”这两段序列截取下来,比较碱基组成是否一致。samtools是个非常不错的工具,可以在这里用于截取fasta文件中特定位置的序列。综上,我们就明确了基因“rpl16”在我们的基因组中的准确位置“join(182615..182623, 183624..184022)”,然后在GeSeq注释结果gbk文件中修改即可。对于蛋白编码的氨基酸序列,如果基因核酸序列组成和参考基因组一致,这种情况下直接复制参考基因组该基因编码的氨基酸序列就可以了(近缘物种叶绿体基因组非常保守,因此绝大多数情况下是一样的);如果有几个碱基不一致,这时还得借助其它工具完成翻译(基因组同源度不高时需要注意,这时候还可能产生提前终止等特殊的情形)。
好了,上述展示了其中一种需要手动确定基因位置的情况,更多情况因篇幅原因就不再举例了。总之记得,GeSeq的原始注释结果一定要检查的,可谓是小问题多多。这个过程还是非常麻烦的,而且需要极大的耐心……最后,你得到了一个校正后的Genbank注释文件,你确定无误了后,就可以将它作个转化了。比如说转成gff文件、fasta文件、tbl文件等,以用于后续特殊的需要。这种情况下,只对的Genbank注释文件作调整,然后再转化得到其他文件还是蛮省事的,这样就不用再单独对GeSeq注释所得的原始gff等文件作调整了。
BioPerl、Biopython等工具包提供了这类的转化脚本,使用起来很方便,就不多提了。大不了使用前搜素一下“genbank转tbl”之类的,网上的方法多的很。
GeSeq的其它功能,如基因组圈图绘制
GeSeq网站中还提供了绘制细胞器基因组圈图的界面,GeSeq绘制的圈图也非常的漂亮。
GeSeq基因组圈图绘制界面:https://chlorobox.mpimp-golm.mpg.de/OGDraw.html
在该界面中,我们上传一个叶绿体/线粒体基因组gbk文件,然后调整好作图参数后,提交,等待一会儿便可得到圈图结果。如下示例。
示例圈图结果如下,还挺好看的吧。如果你也觉得这张图挺漂亮的,那么倒也不必再费时间额外使用circos、circlize、cgview等工具绘制基因组圈图了。
友情链接